#include <iostream>
#include <vector>
#include <fstream>
#include "Import.h"
#include "Node.h"
#include "Material.h"
#include "Property.h"
#include "Quad4.h"
#include "SPC.h"
#include "Force.h"
#include "ImportDisplacements.h"
#include <windows.h>
#include <string>

using namespace std;

int main(int argc, char* argv[]) {

	cout << "#############################################" << endl;
	cout << "####  CPPFEM 1.0                         ####" << endl;
	cout << "####  24.04.14                           ####" << endl;
	cout << "####                                     ####" << endl;
	cout << "####  ACM12                              ####" << endl;
	cout << "####  Stefan Hermann                     ####" << endl;
	cout << "####  Benjamin Scheunert                 ####" << endl;
	cout << "####  Giang To                           ####" << endl;
	cout << "#############################################" << endl << endl;

	Import myImport;

	string sFileName = argv[1];
	cout << "Starting import" << endl;

	myImport.startImport(sFileName);
	///////// TEST   /////////////////////
	vector<Node> nAllNodes;
	vector<Quad4> nAllElements;
	vector<SPC> nAllSPCs;
	vector<Force> nAllForces;

	//Create Material
	Material material1(1, 200000., .3);
	//Create Property
	Property prop1(&material1, 1, 2, 1);

	cout << "Creating nodes" << endl;
	for (int i = 0; i < myImport.iNumberOfNodes; i++) {
		nAllNodes.push_back(Node(i + 1, myImport.grid_X[i], myImport.grid_Y[i], myImport.grid_Z[i]));
	}

	cout << "Creating elements" << endl;
	for (int i = 0; i < myImport.iNumberOfElements; i++) {
		//cout << nAllNodes[myImport.quad4_g1[i]].iId << endl;
		nAllElements.push_back(
				Quad4(myImport.quad4_eid[i], &prop1, &nAllNodes[myImport.quad4_g1[i] - 1], &nAllNodes[myImport.quad4_g2[i] - 1], &nAllNodes[myImport.quad4_g3[i] - 1],
						&nAllNodes[myImport.quad4_g4[i] - 1]));
	}

	cout << "Creating SPCs" << endl;
	for (unsigned int i = 0; i < myImport.spc_displacement.size(); i++) {
		nAllSPCs.push_back(SPC(myImport.spc_sid[i], myImport.spc_node[i], myImport.spc_dof[i], myImport.spc_displacement[i]));
	}

	cout << "Creating forces" << endl;
	for (unsigned int i = 0; i < myImport.force_cid.size(); i++) {
		nAllForces.push_back(Force(myImport.force_sid[i], myImport.force_node[i], myImport.force_cid[i], myImport.force_scale[i], myImport.force_fx[i], myImport.force_fy[i], myImport.force_fz[i]));
	}

	//assign equation numbers
	cout << "Assigning equation numbers" << endl;
	int iCurrentEquationNumber = 0;
	for (unsigned int i = 0; i < nAllNodes.size(); i++) {
		nAllNodes[i].iEquationNr = iCurrentEquationNumber;
		//count number of active dofs
		int iNumberOfDofs = 0;
		for (int j = 0; j < 6; j++) {
			if (nAllNodes[i].aDOFs[j])
				iNumberOfDofs++;
		}
		//cout << "Node number " << nAllNodes[i].iId << " has equation number " << iCurrentEquationNumber << endl;
		iCurrentEquationNumber += iNumberOfDofs;
	}

	cout << "Assembling full stiffness matrix" << endl;
	//Assemble Stiffness Matrix
	vector<unsigned int> vRow;
	vector<unsigned int> vCol;
	vector<double> vCoefficient;

	for (int unsigned iElementCounter = 0; iElementCounter < nAllElements.size(); iElementCounter++) {
		for (int j = 0; j < nAllElements[iElementCounter].iNumberOfDofs; j++) {
			for (int k = j; k < nAllElements[iElementCounter].iNumberOfDofs; k++) {
				//Get equation numbers
				unsigned int iRowEquationNumber = nAllElements[iElementCounter].nNodes[j / 2]->iEquationNr + j % 2; //TODO --> only works for 2 dofs per node
				unsigned int iColEquationNumber = nAllElements[iElementCounter].nNodes[k / 2]->iEquationNr + k % 2; //TODO --> only works for 2 dofs per node
				//Make sure, that it is an upper right coefficient
				if (iRowEquationNumber > iColEquationNumber) {
					unsigned int iTemp = iRowEquationNumber;
					iRowEquationNumber = iColEquationNumber;
					iColEquationNumber = iTemp;
				}

				//cout << "Equation number for row element: " << iRowEquationNumber << endl;

				//Get coefficient
				double dCoefficient = nAllElements[iElementCounter].aStiffnessMartix[j][k];

				bool bFound = false;
				bool bMappingFound = false;

				//Check if coefficient already exists
				for (unsigned int l = 0; l < vCoefficient.size(); l++) {
					if (vRow[l] == iRowEquationNumber) {
						bMappingFound = true;
						if (vCol[l] == iColEquationNumber) {
							vCoefficient[l] += dCoefficient;
							bFound = true;
							break;
						}
					}
				}

				if (!bMappingFound) {
					string sDOF;
					if ((iRowEquationNumber % 2) == 1)
						sDOF = "UY";
					else
						sDOF = "UX";
					//cout << k / 2 << endl;
					//cout << "Mapping: EQ-Number: " << nAllElements[iElementCounter].nNodes[k / 2]->iId << " " << sDOF << endl;
				}

				if (!bFound) {
					vRow.push_back(iRowEquationNumber);
					vCol.push_back(iColEquationNumber);
					vCoefficient.push_back(dCoefficient);
					//cout << "Coefficient " << iRowEquationNumber << "|" << iColEquationNumber << " not found! Will be created." << endl;
				}
			}
		}
	}

	///// Apply DOFs
	cout << "Applying SPCs" << endl;
	for (unsigned int i = 0; i < nAllSPCs.size(); i++) {
		for (unsigned int j = 0; j < nAllNodes.size(); j++) {
			if (nAllSPCs[i].iNodeNumber == nAllNodes[j].iId) {
				for (int k = 1; k <= 2; k++) {
					if (nAllSPCs[i].dofConstrained(k)) {
						nAllNodes[j].aDOFs[k - 1] = false;
						//cout << "Node #" << nAllSPCs[i].iNodeNumber << " DOF " << k << " constrained" << endl;
					}
				}
			}
		}
	}

	//// Write Constrained DOF vector
	vector<unsigned int> vConstrainedDofs;
	//cout << "DEBUG: Constrained DOFs" << endl;
	for (unsigned int i = 0; i < nAllNodes.size(); i++) {
		for (int j = 0; j < 2; j++) {
			if (!nAllNodes[i].aDOFs[j]) {
				vConstrainedDofs.push_back(nAllNodes[i].iEquationNr + j);
				//cout << nAllNodes[i].iEquationNr + j << endl;
			}
		}
	}

	//Reduce Stiffness Matrix
	unsigned int iRedurcedDofs = 0;
	unsigned int iTempEquationNumber = 0;
	for (unsigned int i = 0; i < nAllNodes.size(); i++) {
		for (int j = 0; j < 2; j++) {
			if (!nAllNodes[i].aDOFs[j]) {
				//iTempEquationNumber = nAllNodes[i].iEquationNr + j;

				iRedurcedDofs++;
				for (unsigned int k = 0; k < vCoefficient.size(); k++) {
					if (vRow[k] == iTempEquationNumber || vCol[k] == iTempEquationNumber) {
						//cout << "Erased element " << k << endl;
						vRow.erase(vRow.begin() + k);
						vCol.erase(vCol.begin() + k);
						vCoefficient.erase(vCoefficient.begin() + k);
						k--;
					}
				}
				iTempEquationNumber++;
			} else {
				for (unsigned int k = 0; k < vCoefficient.size(); k++) {
					if (vRow[k] == iTempEquationNumber)
						vRow[k] -= iRedurcedDofs;
					if (vCol[k] == iTempEquationNumber)
						vCol[k] -= iRedurcedDofs;
				}
				iTempEquationNumber++;
			}
		}
	}

	cout << endl << "Writing stiffness matrix, mapping and right hand side to file." << endl;

	///// Write Stiffness Matrix and Mapping
	ofstream fStiffnessMatrix;
	ofstream fMapping;
	ofstream fRightHandSide;

	fStiffnessMatrix.open("stiffnessMatrix.txt", ios::out);
	fMapping.open("mapping.txt", ios::out);
	fRightHandSide.open("rightHandSide.txt", ios::out);

	fStiffnessMatrix.precision(8);
	fRightHandSide.precision(8);

	fStiffnessMatrix << "CPPFEM - Stiffness matrix";

	//Stiffness Matrix
	for (unsigned int i = 0; i < vCoefficient.size(); i++) {
		fStiffnessMatrix << endl << vRow[i] << "\t" << vCol[i] << "\t" << vCoefficient[i];
	}

	//Mapping
	for (unsigned int i = 0; i < nAllNodes.size(); i++) {
		if (nAllNodes[i].aDOFs[0])
			fMapping << nAllNodes[i].iId << "	UX" << endl;
		if (nAllNodes[i].aDOFs[1])
			fMapping << nAllNodes[i].iId << "	UY" << endl;
	}

	//Right Hand Side
	fRightHandSide << "CPPFEM - Right Hand Side";
	bool bFound = false;
	for (unsigned int i = 0; i < nAllNodes.size(); i++) {
		for (int j = 0; j < 2; j++) {
			if (nAllNodes[i].aDOFs[j]) {
				for (unsigned int k = 0; k < nAllForces.size(); k++) {
					if (nAllForces[k].iNode == nAllNodes[i].iId) {
						if (j == 0)
							fRightHandSide << endl << nAllForces[k].fFx;
						else
							fRightHandSide << endl << nAllForces[k].fFy;
						bFound = true;
						break;
					}
				}
				if (!bFound)
					fRightHandSide << endl << "0" ;
				bFound = false;
			}
		}
	}

	fRightHandSide.close();
	fStiffnessMatrix.close();
	fMapping.close();

	cout << "Calling solver" << endl;

	////////////////////////////////////////////////////////////////////////
	STARTUPINFO si;
	PROCESS_INFORMATION pi;

	ZeroMemory( &si, sizeof(si) );
	si.cb = sizeof(si);
	ZeroMemory( &pi, sizeof(pi) );

	// Start the child process.
	if (!CreateProcess( NULL,   // No module name (use command line)
			"Solver_Project.exe",        // Command line
			NULL,           // Process handle not inheritable
			NULL,           // Thread handle not inheritable
			FALSE,          // Set handle inheritance to FALSE
			0,              // No creation flags
			NULL,           // Use parent's environment block
			NULL,           // Use parent's starting directory
			&si,            // Pointer to STARTUPINFO structure
			&pi)           // Pointer to PROCESS_INFORMATION structure
			) {
		printf("CreateProcess failed (%d).\n", GetLastError());
	}

	// Wait until child process exits.
	WaitForSingleObject(pi.hProcess, INFINITE);

	//ShellExecute( NULL, NULL, "calc.exe", NULL, NULL, SW_SHOW );
	//system("Solver_Project.exe");

	/*///// DEBUG Write Stiffness matrix
	for (unsigned int i = 0; i < vCoefficient.size(); i++) {
		cout << vRow[i] << "\t" << vCol[i] << "\t" << vCoefficient[i] << endl;
	}

	cout << endl << endl;
	/// DEBUG Write Mapping
	for (unsigned int i = 0; i < nAllNodes.size(); i++) {
		if (nAllNodes[i].aDOFs[0])
			cout << nAllNodes[i].iId << "	UX" << endl;
		if (nAllNodes[i].aDOFs[1])
			cout << nAllNodes[i].iId << "	UY" << endl;
	}*/

	cout << "Importing displacements" << endl;
	/// Import Displacements
	ImportDisplacements displacementsImport;
	displacementsImport.startImport();

	int iDisplacementsCounter = 0;

	for (unsigned int i = 0; i < nAllNodes.size(); i++) {
		for (int j = 0; j < 2; j++) {
			if (nAllNodes[i].aDOFs[j]) {
				nAllNodes[i].fDisplacements[j] = displacementsImport.vDisplacements[iDisplacementsCounter];
				cout << displacementsImport.vDisplacements[iDisplacementsCounter] << endl;
				iDisplacementsCounter++;
			}
		}
	}

	cout << "Calculating strains and stresses" << endl;
	// Calculate Element Results
	for (unsigned int i = 0; i < nAllElements.size(); i++) {
		nAllElements[i].calculateResults();
	}

	cout << "Writing result file" << endl;
	/// Write results file
	ofstream fResults;

	fResults.open("CPPFEM_result.txt", ios::out);
	fResults << "ALTAIR ASCII FILE" << endl;
	fResults << "$TITLE = Static analysis" << endl;
	fResults << "$SUBCASE = 1	Subcase 1" << endl;
	fResults << "$BINDING = NODE" << endl;
	fResults << "$COLUMN_INFO = ENTITY_ID" << endl;
	fResults << "$RESULT_TYPE = Displacement(v)" << endl;
	fResults << "$TIME = 1.0 sec" << endl;

	for (unsigned int i = 0; i < nAllNodes.size(); i++) {
		fResults << nAllNodes[i].iId << "\t" << nAllNodes[i].fDisplacements[0] << "\t" << nAllNodes[i].fDisplacements[1] << "\t 0" << endl;
	}

	fResults << "$BINDING = ELEMENT" << endl;
	fResults << "$COLUMN_INFO = ENTITY_ID" << endl;
	fResults << "$RESULT_TYPE = Element strains(t), Element stresses(t)" << endl;
	fResults << "$TIME = 1.0 sec" << endl;

	for (unsigned int i = 0; i < nAllElements.size(); i++) {
		fResults << nAllElements[i].iId << "\t" << nAllElements[i].aStrains[0] << "\t" << nAllElements[i].aStrains[1] << "\t" << 0 << "\t" << nAllElements[i].aStrains[2] << "\t" << 0 << "\t" << 0
				<< "\t" << nAllElements[i].aStresses[0] << "\t" << nAllElements[i].aStresses[1] << "\t" << 0 << "\t" << nAllElements[i].aStresses[2] << "\t" << 0 << "\t" << 0 << endl;
	}

	fResults.close();
	return 0;
}
